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FOREWORD 


A  description  of  simple,  dense  fluids  based  on  intermolecular  forces  has 
been  investigated  utilizing  liquid-state  perturbation  theory.  This  work  is 
applicable  to  the  modeling  of  reaction  products  of  explosives.  Funding  for  this 
work  was  provided  by  both  the  Independent  Research  and  Explosives  6.2  Block 
Programs  at  the  Naval  Surface  Weapons  Center. 
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INTRODUCTION 


A  description  of  even  simple  liquids  in  the  high  pressure  domain  is  a 
problem  of  extreme  complexity.  Aside  from  the  work  of  Ross^,  who  employs  a 
variational  technique,  most  efforts  have  been  directed  towards  lattice  theory 
approximations;  but  it  is  well  known  that  the  latter  are  inadequate.  However, 
perturbation  theory  appears  to  provide  a  simple,  but  accurate  characterization 
of  dense  liquids. 

2 

The  formalism  of  Weeks,  Chandler  and  Anderson  is  employed  to  describe 
fluids  whose  molecules  interact  via  a  spherically  symmetric  modified  Buckingham 
(exp-6)  potential.  It  is  assumed  that  the  repulsive  forces  provide  the  dominate 
contribution  to  the  properties  of  dense  fluids.  The  intermolecular  potential 
is  divided  into  a  reference  part  which  is  repulsive  in  character  and  a 
perturbation  which  is  attractive.  Division  of  the  potential  is  carried  out 
so  that  the  resulting  free  energy  is  a  minimum  or  essentially  constant  with 
respect  to  the  break  point  as  suggested  by  Ree3.  With  the  prescription  of 
Verlet  and  Weis  (VW),4  the  radial  distribution  function  for  the  reference  system 
is  akin  to  one  for  hard  spheres  with  a  diameter  which  is  a  function  of  the 
temperature  and  density.  First-order  properties  can  be  written  analytically 
aside  from  several,  simple  numerical  integrals. 

Within  the  above  framework  several  calculations  are  performed.  The 
isothermal  compressibility  and  internal  energy  are  computed  for  a  system 
governed  by  an  exp-6  intermolecular  potential  for  high  temperatures  and 
pressures  and  are  compared  to  Monte  Carlo  (MC)  predictions^.  Then  properties 
for  a  real  system,  liquid  nitrogen,  are  considered.  It  is  assumed  that  the 
internal  modes  are  unchanged  from  the  gaseous  state.  Retention  of  the  free 
rotation  is  based  on  the  observations  of  Smith,  et  al°  which  indicates  that 
there  is  a  transition  from  free  rotation  in  the  solid-e  phase  to  hindered 
rotation  in  the  solid-a  phase. 

*Ross,  M. ,  Phys.  Rev.  A8,  1466  (1973). 

3Weeks,  J.  D. ,  Chandler,  D.  and  Anderson,  H.  C.,  J.  Chem.  Phys.  54_,  5237  (1971). 
3Ree,  F.  H.,  J.  Chem.  Phys.  64,  4601  (1976). 

4Verlet,  L.,  and  Weis,  0.,  Phys.  Rev.  A5.,  939  (1972). 

5Ross,  M.  and  Alder,  B.  J.,  J.  Chem.  Phys.  46,  4203  (1967). 

^Smith,  A.  L.,  Keller,  W.  E.  and  Johnston,  H.  L.,  Phys.  Rev.  7£,  728  (1950). 
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Moderate  temperature  isotherms  and  a  shock  Hugoniot  are  calculated  and 
are  compared  to  experimental  data' *8. 


THEORY 


The  perturbation  theory  to  be  discussed  here  is  based  on  the  assumption 
that  the  repulsive  intermolecular  forces  provide  the  dominant  characteristics  of 
the  material.  This  technique,  originated  by  Weeks,  et  alz,  should  be  quite 
applicable  to  the  high  pressure  domain. 

It  is  assumed  that  the  intermolecular  interaction  is  spherically  symmetric 
and  is  taken  as 


<t>(r)  =  vQ{r)  +  w(r)  ,  (1) 

where  v0(r)  is  the  reference  potential  with  repulsive  character  and  w(r)  is 
the  attractive  perturbation.  The  reference  potential  is  given  by 


vQ(r)  *  - 

4>(r)  -  4>(x) 

rsX 

(2) 

0 

r>\ 

and  the  perturbation  is  taken  as 

( 

r<X 

w(r)  =  j 

i >(r) 

r>x  * 

(3) 

where 

/ 

\  j 

<f>(r)  =  e|(6/a)exp[a(l-r/rm) 

|  -  (yr)6V(l-6/ct).  , 

(4) 

The  potential  parameters  denote  the  usual  quantities.  Specifically,  e 
is  the  well  depth,  a  is  the  steepness  parameter,  and  rm  is  the  position  of  the 
minimum  of  the  potential.  Historically,  the  potential  is  divided  at  the 
minimum2*^;  however,  Ree3  demonstrated  that  the  region  of  applicability  could  be 
extended  to  the  high  density  domain  with  the  additional  flexibility  of  a 
variable  x.  Hence  x  is  here  taken  as  a  variational  parameter. 


7Robertson,  S.  L.  and  Bann,  Jr.,  S.  E. ,  J.  Chem.  Phys.  50_,  4560  (1969). 
8Dick,  R.  D.,  J.  Chem.  Phys.  52,  6021  (1970). 
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q 

Csing  standard  techniques  the  free  energy  can  be  written  as  a  perturbative 
series  in  inverse  temperature,  which  is  well  suited  for  the  shock  wave  domain. 
The  result  for  the  excess  free  energy  per  particle  is 

f  =  f0  +  p/2  /  dr>/(r)g0(r)  +  ...  ,  (5) 

where  f  and  g  (r)  are  the  free  energy  per  particle  and  radial  distribution 
function  for  tRe  reference  system,  respectively,  and  the  number  density  is 
denoted  by  p.  In  this  work  only  the  two  leading  terms  in  the  expansion  of 
Equation  5  are  retained. 

The  steep,  discontinuous  behavior  of  v  (r)  is  reminiscent  of  a  hard  sphere 
potential.  Application  of  a  functional  TayTor  expansion10  yields  the  result 
that  the  free  energy  of  the  reference  system  is  identical  to  second-order 
with  that  for  hard  spheres  with  a  temperature  and  density-dependent  hard 
sphere  diameter,  so  that 


where 

90(r)  =  e"Svo^r^yHs(r)  ,  (7) 

and 

9Hs(r)  *  e*s<^Hs  ^  yHs(r)  • 

The  hard-sphere  diameter  is  determined  from  the  condition 

/d?{e-Svo<r>  -  e-S*HsW  }yHs(r)  -  0  .  (9) 

This  is  the  requirement  that  the  long  wavelength  structure  factors  for  the 
reference  system  and  the  hard-sphere  system  be  equal. 


Combination  of  the.results  of  VW  with  the  Percus-Yevic  description  of 
the  hard-sphere  system11  yields  analytic  results  for  the  thermodynamic  state 
variables  aside  from  several  simple,  numerical  integrals.  The  compressibility 
and  excess  internal  energy  per  particle  are  obtained  from  analytic  derivatives 
of  the  free  energy  as 

z  • »  ■ 1  *  40,  lI0) 

a_^ _ e  ■  . 

9Zqanzig,  R.  W. ,  J.  Chem.  Phys.  22_,  1420  (1954). 

^Anderson,  H.  C. ,  Chandler,  D.  and  Weeks,  J.  D. ,  J.  Chem.  Phys.  56.,  3812  (1972). 
nWertheim,  M.  S.,  Phys.  Rev.  Lett.  10,  E501  (1963). 


6 


NSWC  TR  80-419 


It  was  found  that  the  compressibility  defined  in  this  manner  rather  than 
the  original  analysis  of  VW  compared  more  favorably  with  MC  results  in  higher 
pressure  and  temperature  region. 


RESULTS  AND  DISCUSSION 


The  isothermal  properties  for  a  typical  exp-6  potential  with  the  parameters 
ct  =  13.5  i 
rm  =  3.85  x  10* 10  m  , 
e/k  =  122°K  , 

are  given  in  Table  1.  Here,  the  reduced  temperature  is  defined  as  T*  =  kT/e,  and 

the  reduced  density  is  taken  as  p*  =  (r^yf2)p.  For  a  material  such  as  argon,  the 

lower  temperature  isotherm  ranges  from  2.5  GPa  to  16.6  GPa,  while  the  higher 
temperature  isotherm  varies  from  16  GPa  to  100  GPa.  It  is  seen  that  perturbation 
theory  is  in  excellent  agreement  with  Monte  Carlo  calculations,  except  for  the 
very  high  density  regions  for  the  respective  isotherms.  In  both  instances  these 
larger  deviations  result  from  the  entrance  into  the  meta  stable  region  for  the 
hard-sphere  reference  system. 

To  consider  a  real  system,  moderate  temperature  isotherms  for  liquid 
nitrogen  are  computed.  The  values  of  the  potential  parameters  are  taken  as 

o  =  13.6  , 

rm  =4.12  x  10'10m  , 

e/k  =  100. 0°K 

The  compressibility  as  a  function  of  the  pressure  is  shown  in  Figure  1. 

The  solid  curves  are  the  perturbative  predictions,  and  the  circles  represent 
the  experimental  data  It  is  seen  that  the  calculation  exhibits  very  good 
agreement  with  the  data,  with  a  maximum  error  of  5".  Thus,  the  theory 
describes  the  liquid  quite  well  in  the  moderate  pressure  and  temperature  range. 

As  a  final  application  of  this  technique,  the  shock  compression  of 
liquid  nitrogen  is  considered.  The  initial  conditions  for  the  material  are 
given  as  T0=77°K  and  p0=8.1  x  102  kg/m2.  The  same  potential  parameters  used 
in  the  previous  calculation  of  moderate  temperature  isotherms  are  used  for 
this  computation.  In  this  calculation  a  value  of  volume  for  the  compressed 
state  was  chosen,  and  then  a  value  of  the  temperature  was  found  by  iteration 
to  satisfy  the  energy  conservation  condition, 

E  -  Eq  =  (P+P0)(v0-v)  ,  (12) 
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TABLE  1 

THERMODYNAMIC  PROPERTIES  FOR  EXP. -6  POTENTIAL 


T*  =  20 


Monte  Carlo5 

Perturbation  Theory 

p* 

Z 

BE 

Z 

be 

.9 

3.27 

.28 

3.30 

.28 

1.0 

3.81 

.38 

3.83 

.38 

1.25 

5.52 

.73 

5.56 

.75 

1.50 

8.07 

1.34 

7.85 

1.32 

1.75 

11.34 

2.17 

10.3 

2.13 

T*  =  100 

1.264 

3.06 

.55 

3.09 

.57 

1.473 

3.71 

.74 

3.71 

.74 

1.768 

4.81 

1.08 

4.75 

1.06 

2.431 

8.17 

2.20 

7.57 

2.03 

2.701 

9.86 

2.78 

9.11 

2.78 
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0.2  0.4  0.6  0.8  1.0 

P(GPa) 


FIGURE  1  COMPRESSIBILITY  VS.  PRESSURE  FOR  LIQUID  NITROGEN 
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where  the  zero  subscript  refers  to  material  in  the  undisturbed  state  ahead 
of  the  shock  front.  In  Equation  (12),  E  and  v  are  the  total  specific  energy 
and  specific  volume,  respectively. 

The  theoretical  shock  Hugoniot  in  the  pressure-volume  plane  for  liquid 
nitrogen  is  compared  with  the  experimental  data®  in  Figure  2.  The  data 
could  be  interpreted  as  having  a  cusp  in  the  14  GPa  region,  signalling 
the  onset  of  a  phase  transformation;  however,  there  is  good  agreement 
between  the  present  calculation  and  experiment.  Similar  results  are 
exhibited  in  Figure  3  in  the  shock  velocity  particle  velocity  plane.  The 
decline  in  the  slope  is  well  characterized  by  the  single-phase  computation. 
For  completeness,  the  theoretical  Hugoniot  is  presented  in  Table  2. 


SUMMARY 

Thermodynamic  properties  for  high  temperatures  and  pressures  have  been 
calculated  using  an  exp-6  potential  and  liquid-state  perturbation  theory. 

The  results  were  in  good  agreement  with  MC  computations,  except  for  extremely 
dense  states  which  are  close  to  the  phase  line  of  the  hard  sphere  reference 
system.  Theoretical  predictions  for  the  compressibil ity  of  liquid  nitrogen 
for  moderate  temperatures  were  in  accord  with  the  experimental  data.  The 
calculated  shock  Hugoniot  for  liquid  nitrogen  was  also  in  close  agreement 
with  experimental  results.  It  is  apparent  that  the  theoretical  technique 
utilized  in  the  calculations  provides  a  convenient  and  accurate  description 
of  simple  fluids  in  the  moderate  as  well  as  high-pressure  domains. 


10 


V/V0 

FIGURE  2  SHOCK  PRESSURE  VS.  RELATIVE  VOLUME  FOR  LIQUID  NITROGEN 


NSWC  TR  80-419 


TABLE  2 

SHOCK  HUGONIOT  FOR  LIQUID  NITROGEN 


v/vo 

P(GPa) 

.653 

1.86 

.628 

2.47 

.570 

4.83 

.531 

7.56 

.501 

10.7 

.494 

11.6 

.468 

16.0 

.448 

20.5 

.432 

25.3 

.419 

30.0 

.409 

34.2 

.400 

39.0 

T(°K) 

U(km/ 

316.2 

2.56 

418.3 

2.84 

874.7 

3.70 

1446 

4.44 

2137 

5.12 

2348 

5.30 

3364 

6.05 

4486 

6.73 

5714 

7.37 

6977 

7.94 

8134 

8.41 

9467 

8.90 

)  Up (km/sec) 

0.89 

1.06 

1.59 

2.08 

2.55 

2.68 

3.22 

3.72 

4.19 

4.61 

4.97 

5.34 
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OAK  TECHNICAL  REPORTS  PLEASE  COMPLETE  THE 
FORM  BELOW: 

TO  ALL  HOLDERS  OF  NSWC/TB  80-419 

by  Hermenzo  D.  Jones,  Code  R13 

DO  NOT  RETURN  THIS  FORM  IF  ALL  INFORMATION  IS  CURRENT 

A.  FACILITY  NAMC  ANO  ADDRESS  (OLD)  (Show  Zip  Code) 


NEW  AOORESS  (Show  Zip  Code) 


■.  ATTENTION  LINE  ADDRESSES: 


C. 


|  I  REMOVE  THIS  FACILITY  FROM  THE  DISTRIBUTION  L.ST  FOR  TECHNICAL  REPORTS  ON  THIS  SUBJECT. 


D. 

number  of  copies  desired 
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